B. rapa Gene Networks

Julin N Maloof

July 17, 2018

Intro

Background

  • Goal: find genes & networks connected to B. rapa growth
  • Rob Baker, Cynthia Weinig, Steve Welch, etc collected growth data on B. rapa RIL set in field setting.
  • Growth functions were fitted for each line.
  • We have gene expression for the same lines.
  • Are there genes whose expression is associated with growth?
  • Can they help explain growth variation between lines?

Approach

  • Use Mutual Rank co-expression networks:
    • For each growth trait, build a network of genes that have a similar pattern across the RILs
    • In this case, look for association by mutual rank (MR) (next slide)
  • Alternatively use WGCNA:
    • Build WGCNA co-expression network independent of growth traits
    • Then look for correlation between network module Eigen genes and growth traits

Mutual Rank

Basic idea of Mutual Rank

  1. Build a correlation matrix of gene expression and trait values
  2. For each gene or trait, rank its correlations with the other genes and traits. The gene with the highest correlation gets rank 1, etc.
  3. For each gene x gene, gene x trait, or trait x trait pair, calculate the geometric average. This is the mutual rank.
  4. Look at the network of genes (and traits) around each trait.
  5. Determine significance by permutation.

Starting data

Line UNr12 UNHmax12 UNiD12 UNd12 UNSTP12 Bra016271 Bra038030 Bra040132 Bra017243
RIL_1 0.0127 40.1 659.2 944.4 38.4 8.8 5.7 4.8 4.7
RIL_103 0.0141 47.0 617.6 870.7 45.6 8.9 1.2 4.8 -3.0
RIL_115 0.0135 26.3 652.6 910.2 24.7 8.9 0.9 -1.6 -3.2
RIL_12 0.0145 19.9 516.6 762.6 18.8 2.4 -1.6 5.0 -3.9
RIL_124 0.0120 44.7 695.4 966.1 41.8 9.2 5.9 -1.1 5.0
RIL_131 0.0138 28.0 487.3 750.1 26.9 2.9 5.9 4.9 -3.3
RIL_136 0.0107 34.0 575.6 884.7 32.2 8.9 6.1 -1.6 5.4
RIL_143 0.0128 39.1 742.0 997.9 36.6 2.6 5.8 -2.4 5.3
RIL_146 0.0113 37.8 596.4 896.6 35.8 2.6 5.6 -2.1 4.7
RIL_147 0.0140 38.2 693.8 936.8 36.0 2.2 -0.7 -1.7 5.2

Build correlation matrix

UNr12 UNHmax12 UNiD12 UNd12 UNSTP12 Bra016271 Bra038030 Bra040132 Bra017243
UNr12 1.00 -0.25 -0.32 -0.51 -0.23 -0.04 -0.39 0.08 -0.15
UNHmax12 -0.25 1.00 0.79 0.77 1.00 0.09 0.28 0.13 -0.01
UNiD12 -0.32 0.79 1.00 0.97 0.78 -0.03 0.30 -0.09 0.20
UNd12 -0.51 0.77 0.97 1.00 0.75 -0.03 0.36 -0.10 0.21
UNSTP12 -0.23 1.00 0.78 0.75 1.00 0.09 0.27 0.14 -0.02
Bra016271 -0.04 0.09 -0.03 -0.03 0.09 1.00 -0.01 -0.02 -0.03
Bra038030 -0.39 0.28 0.30 0.36 0.27 -0.01 1.00 0.03 0.04
Bra040132 0.08 0.13 -0.09 -0.10 0.14 -0.02 0.03 1.00 -0.20
Bra017243 -0.15 -0.01 0.20 0.21 -0.02 -0.03 0.04 -0.20 1.00

Calculate mutual ranks

UNr12 UNHmax12 UNiD12 UNd12 UNSTP12 Bra030612 Bra010212 Bra102212 Bra101312 Bra007112
UNr12 1 1694 768 22 2041 6291 9173 3227 2416 9701
UNHmax12 1694 1 3 3 2 9702 6554 1623 8654 3079
UNiD12 768 3 1 2 3 8770 8021 2076 7804 1536
UNd12 22 3 2 1 4 9449 8108 1767 6085 2078
UNSTP12 2041 2 3 4 1 9623 6041 1588 8241 3145
Bra030612 6291 9702 8770 9449 9623 1 2194 7135 7878 9667
Bra010212 9173 6554 8021 8108 6041 2194 1 4916 9702 8193
Bra102212 3227 1623 2076 1767 1588 7135 4916 1 6312 4154
Bra101312 2416 8654 7804 6085 8241 7878 9702 6312 1 5763
Bra007112 9701 3079 1536 2078 3145 9667 8193 4154 5763 1

What do we get?

MR < 20

MR < 30

MR < 50

Significance thresholds

  • Randomize traits relative to gene expression
  • Make a graph, count # of connections
  • Repeat 100 times
  • Calculate 95th percentile # of connections
  • You would expect fewer than this number 95% of the time in random data

Significance thresholds

Permutation analysis: UN network
threshold 0.95 expectation observed
blup.mr.graph.UN.2012.10 10 0 1
blup.mr.graph.UN.2012.20 20 1 13
blup.mr.graph.UN.2012.30 30 4 28
blup.mr.graph.UN.2012.50 50 12 57

So the networks are enriched for “real” connections

Candidates from UN condition

  • ATMES10: Methyl jasmonate esterase activity; increases JA. Negatively correlated with growth and final height.
  • AtAAS: Induced by JA. Positively correlated with growth (strange…)
  • ENH1: Enhancer of Salt Over Sensitive. Negatively correlated.
  • BEL1: interacts with STM and KNAT1. Required for ovule specification, but could be more general growth regulator. Restricts WUS expression. Negatively correlated with growth.
  • LPPEPSILON 2: lipid phosphatidic acid phosphatase. In chlroplast
  • LBD37: LOB domain containing protein. Nitrogen Metabolism. Negatively correlated with growth.

MR genes and growth QTL

Do any of the the MR genes co-localize with QTL peaks for growth?

Yes!

MR Cutoff Network Genes Network Genes under QTL P-Value for enrichment
10 1 0 1
20 13 9 6.9 E-13
30 28 9 5.0 E-09
50 57 21 9.9 E-21

Caveat: Linkage…

Could any of these be causal?

  • Null hypothesis: Causal QTL and expression QTL for MR genes are genetically linked by not causally linked
  • Alternative: Causal QTL and expression QTL are the same gene.

\[ growth \sim QTL\_marker\tag{1} \\ \\ VS \] \[ \\growth \sim QTL\_marker + gene\_expression\tag{2} \]

  • Model 2 is favored (Model1 AIC = -1702; Model2 AIC = -1814)

MR eQTL growth QTL overlap

Do eQTL for the MR genes overlap with the growth QTL?

Any trans eQTL?

WGCNA

Basic idea

  • WGCNA is an alternative method for reconstructing genetic networks
  • Also based on gene co-expression
  • Use a ‘soft threshold’ to build a network that approximates scale-free topology
  • For my network I did NOT include growth trait info.

Correlation with traits

After reconstructing the network, calculate the Eigen gene for each module and then correlate with traits:

GO enrichment of correlated modules

Blue: Cell division (pos correlation)

module term ontology FDR
blue cell division BP 0.0000382
blue nucleosome assembly BP 0.0001382
blue microtubule-based movement BP 0.0010579
blue auxin metabolic process BP 0.0027352
blue regulation of meristem structural organization BP 0.0027352
blue leaf development BP 0.0035562
blue xylem and phloem pattern formation BP 0.0181424
blue microtubule nucleation BP 0.0188533
blue positive regulation of microtubule polymerization BP 0.0188533
blue regulation of cell cycle BP 0.0192759
blue DNA replication initiation BP 0.0229370
blue abaxial cell fate specification BP 0.0229370
blue fruit development BP 0.0229370

GO enrichment of correlated modules

Cyan: Translation (pos correlation)

module term ontology FDR
cyan translation BP 0.0000000
cyan cytosolic large ribosomal subunit CC 0.0000000
cyan nucleolus CC 0.0000000
cyan cytosolic small ribosomal subunit CC 0.0000000
cyan ribosome biogenesis BP 0.0000000
cyan vacuolar membrane CC 0.0000006
cyan mitochondrial respiratory chain complex I CC 0.0000310
cyan cell wall CC 0.0021796
cyan cytosolic ribosome CC 0.0046775
cyan chloroplast CC 0.0261363
cyan plasma membrane CC 0.0349344

GO enrichment of correlated modules

MidnightBlue: Stress/defense? (pos correlation)

module term ontology FDR
midnightblue vacuole CC 0.0000209
midnightblue toxin catabolic process BP 0.0000209
midnightblue response to wounding BP 0.0000499
midnightblue ER body CC 0.0000765
midnightblue sulfate assimilation BP 0.0000974
midnightblue glucosinolate biosynthetic process BP 0.0015902
midnightblue regulation of response to red or far red light BP 0.0042999
midnightblue chloroplast stroma CC 0.0044724
midnightblue male gamete generation BP 0.0044724
midnightblue response to cadmium ion BP 0.0081931
midnightblue post-embryonic root development BP 0.0081931
midnightblue glucosinolate catabolic process BP 0.0081931
midnightblue response to jasmonic acid BP 0.0098597
midnightblue response to herbivore BP 0.0102517
midnightblue regulation of glucosinolate biosynthetic process BP 0.0127414
midnightblue glutathione metabolic process BP 0.0214903
midnightblue cellular response to water deprivation BP 0.0236303
midnightblue heterotetrameric ADPG pyrophosphorylase complex CC 0.0260464
midnightblue stamen formation BP 0.0260464
midnightblue response to misfolded protein BP 0.0260464
midnightblue proteasome core complex assembly BP 0.0260464
midnightblue mitochondrial proton-transporting ATP synthase complex CC 0.0260464
midnightblue response to salt stress BP 0.0260464
midnightblue nitrile biosynthetic process BP 0.0260904
midnightblue oxidation-reduction process BP 0.0346123
midnightblue peroxisome CC 0.0346123
midnightblue defense response to bacterium BP 0.0346123

GO enrichment of correlated modules

Brown: ? (neg correlation)

module term ontology FDR
brown protein serine/threonine phosphatase complex CC 0.0005793
brown actin cytoskeleton organization BP 0.0005793
brown protein dephosphorylation BP 0.0019447
brown small GTPase mediated signal transduction BP 0.0495296

Module Eigen gene eQTL

eQTL / trait QTL overlap

A10 (cyan, midnightblue, blue) and A03 (darkslateblue) overlap with growth QTL

Summary

Summary

Network analysis has defined possible mechanisms underlying growth trait QTL.

Linkage is a major caveat

Model comparison is encouraging that more than linkage is going on